Mon. Not. R. Astron. Soc. 000, 000-000 (0000) 



Printed 1 February 2008 



(MN style file v2.2) 



Production of trans-Neptunian binaries through 
chaos-assisted capture 

Ernestine A. Lee/* Sergey A. Astakhov,^'^ f and David Farrelly^J 

^FivePrime Therapeutics, 1650 Owens Street, Suite 200, San Francisco, CA 94158-2216, USA 
^Department of Chemistry and Biochemistry, Utah State University, Logan, UT 84322-0300, USA 
^ UniquelCs, Stroiteley 1, Saratov, 4 1 0044i Russia^ 



ABSTRACT 

The recent discovery of binary objects in the Kuiper-belt opens an invaluable 
window into past and present conditions in the trans-Neptunian part of the 
Solar System. For example, knowledge of how these objects formed can be used 
to impose constraints on planetary formation theories. We have recently pro- 
posed a binary-object formation model based on the notion of chaos-assisted 
capture. In this model two potential binary partners may become trapped for 
long times inside chaotic layers within their mutual Hill sphere. The binary 
is then captured permanently through gravitational scattering with a third 
"intruder" body. The creation of binaries having similarly sized partners is an 
ab initio prediction of the model which also predicts large binary semi-major 
axes and moderately eccentric mutual orbits similar to those observed. Here 
we present a more detailed analysis with calculations performed in the spatial 
(three-dimensional) three- and four-body Hill approximations. It is assumed 
that the potential binary partners are initially following heliocentric Keplerian 
orbits and that their relative motion becomes perturbed as these objects un- 
dergo close encounters. First, the mass, velocity, and orbital element distribu- 
tions which favour binary formation are identified in the circular and elliptical 
Hill limits. We then consider intruder scattering in the circular Hill four-body 
problem and find that the chaos-assisted capture mechanism is consistent with 
observed, apparently randomly distributed, binary mutual orbit inclinations. 
It also predicts asymmetric distributions of retrograde versus prograde orbits. 
The time-delay induced by chaos on particle transport through the Hill sphere 
is analogous to the formation of a resonance in a chemical reaction. Implica- 
tions for binary formation rates arc considered and the "fine-tuning" problem 
recently identified by Noll et al. ( 20071) is also addressed. 

Key words: celestial mechanics - methods: N-body simulations - minor plan- 
ets, asteroids - Kuiper Belt - binaries 



1 INTRODUCTION 

The discovery of the first trans-Neptunian - or Kuiper- 
- binary came in two parts punctuated by a 
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prolonged hiatus; the erstwhile planet Pluto was dis- 
covered by Clyde TombauRh in 1930 (jSliphed [l930l : 
iBasri fc Brown l2006l : ISoteij |2006| ) whereas its compan- 
ion, Charon, remained undetect ed until 1978 when it 
was d iscovered almost by chance (|Christv fc HarringtonI 
Il978l ). The identification of the second binary trans- 
Neptunian object (TNO) was also a two-step - and sim- 
ilarly fortuitous - process; during the analysis of orbit 
recovery data for 1998 WW31 it was recognized tha t 
this object is actually a binary (|Veillet et alj I2OO2I ). 
This discovery sparked observing campaigns designed 
to discover other tran s-Neptunian binaries (TNBs) 
(IBrown fc TruiiUol 1200^ and, since then, the inven- 
tory of TNBs has e xpanded rapidly; including Centaurs 
(|No1I et alJ 2006bl). currently about 30 such binaries 



are known (|johnstorJ|2cio7l: INo11||2003I: lOsip et al.ll2003| : 



Farinella et al. 



Kern fc Elliot' '2006b"a'; 'Noll' '200^; 'Noll et al.' '2006a'! 
2007; Stansbcrry ct al. 2006; Stephens & NoU 2006) 

In pl anetary ph ysics, as in other areas of as- 
tronomy Hdl l2003l ). the motivation for finding bi- 
naries is that they can be recruited as natural 
probes of past and present conditions in the ir lo- 
cale - in this in stance the Kuiper-be lt (ISter: 
WiUiamsl Il997|: iKeiivon fc Luul Il998l: iJewit 



199e 



199c 



2OO0I: Malhotra et al.li20qol: lAUen et al 



2001: iKenvonI I2OO2I: ILuu fc JewittI |2002|: ISchulz 
2OO2I : iLevison fc Morbid elli' '20031: IJones et al.' '20051: 
Levi so n et all 



Chiang et all l2006l : ILevi so n et all 120061 : iMorbideUii 
2OO6I I. For example, knowledge of binary mutual or- 
bits provides for a direct and accurate measure- 
ment of the total mass of the system. This then 
paves the way for the determination of other pr oper- 
ties, e.g., bulk densities and mass distribu t ions (jToth 



20061 : 



1999; Noll 'l2003l: lllestroffer fc Vachierll2005l: iNoll et al, 
2006al: iKenvon fc BromlevI |2004 ICruikshank et al 



20061 : iDescamps et al.ll2007l ) . These data, together with 



information about binary frequencies, are key to under- 
standing the format i on and evolution of accr etion disks 

(|Luu fc Jewittll2002l : [Cruikshank et al.ll2006l). 

Amongst Solar System binaries ( Merline et al.l 

I2OO2I : Noll 2003), TNBs are of particular interest be- 
cause of their rather unusual orbital and physical prop- 
erties. These include large, moderately eccentric, mutual 
orbits; randomly distributed inclinations; and a seeming 
preference for the binary partners to have comparable 
sizes (Margot 2002; Bu rns 2004; Noll 2003; Jones_et al. 
I2OO5I : iNoll et al.l [20()6al . [jpOTi ). The discovery of pecu- 
liar properties is a felicitous event because these obser- 



Edgcworth-K uiper- belt, was apparently first sp eculated on 
\w Lconar^ l ll93p[). After late r conjecture s by [ Edeeworthl 
1 I194I il949f) and iKuipe j lll95ll '). iFernandej (I198ot" provided 
the first quantitative discussion of the possibility of a comet- 
belt b e yond Neptune using a Mo n te Car l o met hod - see 
iJewitj ll999h : lLuu fc Jewittj ll2002h : I JewittI ll2007f) . 



vations can be used to tighten constraints on theories 
of the formation and evolution of the Kuiper-belt and, 
more generally, the Solar System. For example, the dis- 
tribution of binary mutual orbit inclinations potentially 
provides insight into the velocity dispersion of bod- 
ies in the primordial Kuiper-belt (IChiang et al] I2OO6I : 
INoU et al.ll2007l : IColdreich et aLllioQ/J ). 

The implication of the discovery of large binary 
orbits with roughly equal mass partners is that TNBs 
did not form through physical collisions; generally these 
are expected to produce objects with ra ther asymmet- 
ric masses and relatively sma l l orbits |Margotl |2002| : 
iBurnsI 12004 iDurda et all |2004| : ISternI I2OO2I ). However, 
an important caveat applies; TNO binaries are diffi- 
cult to observe eve n with the largest ground-based tele- 
scopes (|Tothll 19991 ) which opens up the possibility that 
the apparent preference for large, symmetric-mass bi- 
naries is, in reality, the result of observational bias 
(iBurnslliooi ). Fortunately, it has been possible signif- 
icantly to better characterize these objects with the 
Hubble Sp ace Telescope (HST ) High Resolution Cam- 
era (HRC): [No11 et al.l (|2006al ') have recently made very 
deep observations of TNBs using the HRC in which they 
determined magnitude differences, Amag, between bi- 
nary components. These are the first observations capa- 
ble of measuring the relative frequency of symmetric 
and asymmetric binaries, albeit with the a ssumption 
that re lative brightness is a proxy for size. iNoU et al.l 
(|2006al ) observed statistically significant clustering of bi- 
naries with Amag < 1 and concluded that the preference 
for symmetric binaries is probably real and peculiar to 
TNBs. 

We have recent ly proposed a dynam ical TNB for- 
mation mechanism (|Astakhov et al.ll2005l ) based on the 
idea of ch aos-assisted capture (CAC) in the circular Hill 
prob lem (|Astakhov et all l2003l : lAstakhov fc Farrellvl 
I2OO4I ). Chaos-assisted capture happens because the in- 
terface between bound and scattering regions of phase 
space in the Hill problem consists of "sticky" chaotic lay- 
ers in which trajectories can become trapped and mimic 
regular (i.e., non-chaotic) orbits for lo ng peri ods of time 
(IPerry fc Wiggin s" 1994: S imo fc StuchfcoOOi : IZaslavskvl 
Il985l : lAstakhov et al,,.20o3 r The extension of the life- 
time of the transient binary through its entanglement 
in chaotic layers then provides the opportunity for per- 
manent capture. We proposed that the binary is stabi- 
lized through gravitational scattering with a smaller "in- 
truder" particle. Subsequent intruder scattering events 
gradually reduce the size of the binary orbit and this 
process eventually results in an essentially Keplerian 
mutual orbit. 

Numerical simulations in the CAC model 
(|Astakhov et al.l |2005| ) indicated that symmetric 
binaries (i.e., binaries consisting of similar sized part- 
ners) appear to be created preferentially. This was 
explained as being the result of chaos preferentially 
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destabilizing asymmetric mass binaries as compared to 
symmetric binaries during encounters with intruders. 

Alternative TNB formation models include: phys- 
ical collisions of two objects which then fuse into 
a single object; because all of this is assumed to 
take place inside the Hill spherjfl of a third ob ject a 
binary eventually results JWeidenschillind I2OO2I ): dy- 
namical fri ction ( Goldreich et alj |2002| ): gravitational 
scatter ing jGoldreich et alj 2002^); and exchange reac- 
tions l|Funato et al. 2004 ). D iscus si on of t hese mod- 
els can be found in e.g. iNollI (|2006|): INoII et al.l 
J2OO7'); 'Astakhov et al. ' (I2OO5I ): iKern fc EUiotI (|2006bl ): 
[Cruikshank ct al. (2006). Generally, physical collisions 
alone are unlikely to have formed binaries with large mu- 
tual semi-major axes, in part, because TNBs have signif- 
icantly more an gular momentum than typically results 
from a collision ('Margot"2002':'Burns"2004:'Durd a et al.l 
|2004; Stern 2002; Caimp 2005; Chiang ct al. 20di) al- 
though the satellites of Pluto and of some other TNOs 
likely have a collisional origin (|Canup|2005l : lBrown et al.l 
l2007l : lMorbidellill2007l ). This suggests that perhaps the 
majority of TNBs have a dynamical origin, e.g., in- 
volving gravitational scattering or dynamical friction 
(|Kern fc EUiotI l2006bl ). The CAC model invokes scat- 
tering inside the Hill sphere of three small bodies - the 
potential binary partners (the "primaries") and a third, 
intruder particle. Thus, the overall process is four-body 
including the Sun. 

Here we present a more detailed investigation of the 
initial stages of capture in the CAC scenario, initially 
in the spatial three-body circular and elliptical Hill ap- 
proximations. In part this is an attempt to understand 
how the asymptotic pre-encounter orbital elements de- 
termine capture probabilities. Subsequently, we focus on 
how the mechanism in the four-body circular Hill prob- 
lem, in particular, depends on the masses and velocities 
of the three small bodies. 

A note on terminology is in order: prior to perma- 
nent capture a temporary binary must form. This object 
will be referred to as a transient, quasi-bound, or proto-, 
binary. Immediately after capture the binary typically 
has a very large non-Keplerian orbit; henceforth we re- 
fer to it as a "nascent" binary. The mechanics of or- 
bit reduction of the nascent binary by further intruder 
scattering events is not considered in detail in this arti- 
cle and will be reported separately. For convenience the 
main symbols and terms used and their definitions are 
collected together in Table 1. 

The paper is organized as follows: Section [5] in- 
troduces the Hamiltonian and equations of motion 
of the elliptic spatial Hill problem. Also in Sec. [2] 
we briefly review the CAC mechanism and define 



^ a region wherein mutual gravity dominates solar differen- 
tial gravity - see Table 1 



orbital elements suitable for describing the motion 
of the primaries at infinity in the Hill approxima- 
tion. A similar approach has been employed in stud- 
ies of collisions between planetesimals in a disk re- 



Greenzweig & Lissauerl 


1990[ 1992: WetheriU & StewartI 


1993 : Nakazawa et al. 


1989: Ohtsuki & Ida 1990l: 


Ohtsuki et al.ll2002l:|Stewart & Ida|2000l) or dust grains 


in Dlanetarv rings JPetit & HenonI 1987allbh. Distribu- 



tions of orbital elements at infinity which can lead to 
capture in the spatial three-body circular and elliptical 
Hill approximations are then computed. Four-body in- 
truder scattering simulations in the circular Hill approx- 
imation are described in Sec. |3] Results are presented 
in Sec. [l] Comparison is mad e with the predictions a nd 
assumptions of the models of IColdreich et ahl l|2002h in 
Sec.O The binary hardening mechanism is briefly con- 
sidered in Sec. [6] limitations of our calculations are con- 
sidered in Sec. [7] and conclusions are in Sec. [S] 



2 THREE-BODY HILL APPROXIMATION 

The general plan of attack is as follows: two individual 
TNOs are initially assumed to be orbiting the Sun on (in 
general, elliptical) Keplerian orbits proximate to a cen- 
tral Keplerian orbit with semi-major axis a©. These ob- 
jects may, depending on their relative orbital elements, 
approach to within their mutual Hill sphere radius, Rh 
- see Table 1 (jMurrav fc Dermott|[l999l : iGoldreich et al.l 
I2OO4I ). If the two bodies then happen to get caught up in 
a chaotic layer the resulting quasi-bound binary may be 
permanently captured by gravitational scattering with 
an intruder. 

The circular restricted three-body problem 
(CRTBP) and the three-body HiU problem s (|Hillll 18781 : 
ISzebehelvllT967l : iMurrav fc Dermottlll999l ) have proved 
to be fruitful as test-beds for the study of capture and 
chaos in dynamical s y stems ifMurrav fc Dermot3 19991: 
ISimo & Stuchi' '2000^; 'Belbruno' '2004': 'Asta khov et al.l 
3003; Ast akhov fc Fa rrcUv 2004; Xia 1992). In its most 
usual form, the Hill problem consists of two small 
bodies, mi and m2, orbiting a third, much larger, body, 
mo (hereafter, the Sun, i.e., mo = Mq) with the center 
of mass - "the barycenter" - of mi and m2 following a 
circular orbit. If the barycenter follows an elliptical or- 
bit then the elliptical Hill problem results. Throughout 
we will refer to the candidate binary partners, mi and 
m2, as the "primaries" and will assume, without loss 
of generality, that m2 ^5 mi . Later a fourth body - the 
"intruder," mass ma - will be introduced which may 
undergo gravitational scattering with the binary under 
the governing influence of the Sun. 

In Hill's problem at the Hill sphere radius an equi- 
librium exists between the solar tidal force and the mu- 
tual attraction of the two particles. Thus, the Hill sphere 
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Table 1. Main symbols and terms used and their definitions 



Symbol/Term 



Definition 



G 

Mq 
mi , m2 
m,3 

d ~ 1 g cm~^ 
a0 ~ 45 AU 

60 

/© 



mi+m2 \ 3 
. 3Mq j 

(a, e, i, T, fii, oj) 

6 = |a| = \a2 — ai| 

S 

V 

V 

^ |^ G(rTVi+m2) j 2 

VH ~ 1.2 

V/f ~ 4.4 km/s 

r 

^0 

T = 2iT 

Transient, proto- or quasi-bound binary 

Nascent binary 

Binary hardening (softening) 



Gravitational constant 
Solar mass 

Binary partner masses - the "primaries" 

Mass of fourth body - the intruder - scattered by primaries 

Diameters of primaries 

Physical density of bodies 

Barycenter heliocentric semi-major axis 

Barycenter heliocentric eccentricity 

Barycenter heliocentric true anomaly 

Barycenter heliocentric orbital frequency 



■T^rjz ''^ units 



Radius of binary Hill sphere: Rjj 
Semi-major axes of primaries 
Hill orbital elements of binary barycenter 
Impact parameter 
Surface mass density of primaries 
Velocity dispersion of primaries in physical units 
Velocity dispersion of primaries in Hill units 

Hill velocity of primaries in physical units 
Hill velocity of primaries in Hill units 

Keplerian velocity of primaries at 45 AU in physical units 
Keplerian velocity of primaries in Hill units 
Jacobi constant 

Orbital period at 45 AU in physical units (fa 300 years) 
Orbital period at 45 AU in Hill units 

Hill lifetime; maximum lifetime of an orbit inside the Hill sphere (years) 
Temporary binary inhabiting the Hill sphere 

Newly stabilized binary following a single intruder scattering event 
A process which increases (decreases) the binary binding energy 



radius provides a natural distance scale to describe the 
motion of particles for which solar tides are a pertur- 
bation. Therefore, throughout we use Hill units which 
"blow-up" distances - and orbital elements - in the vicin- 
ity of the binary barycenter so that the radius of the 
Hill sphere, together with distances and velocities a re 
all typically of order unity ijMurrav fc Dermott|[l999l ). 

Although the circular Hill prob lem can be de- 
rived as a special ca se of the CRTBP (|Szebehelvlll"967l : 
iMurrav k, Dermottlll999l), as has been pointed out by 
Heno n and Petit (flKnon fc Petitlll986l : (Petit fc Henonl 
Il986l ). the CRTBP and Hill problems are physically dis- 
tinct. In the CRTBP two masses are assumed to be 
much larger than the third; the Hill problem emerges 
if it is further assumed that mo >> mi >> m2 - this 
"hierarchical" case is shared by the CRTBP and Hill's 



problem. However, Hill's approximation is more gen- 
eral and is valid for arbitrary values of the mass ratio 
mi/m2. Therefore, the Hill approximation is suitable 
for treating the dynamics of TNBs for which the ra- 
tio m\/m2 is, in principle, arbitrary and, remarkably, is 
often of order unity. 



2.1 Hamiltonian and equations of motion 

In practice the three-dimensional (spatial) elliptic 
Hill problem can be derived most directly from 
the elliptic restricted three-body problem in a sim- 
ilar procedure to the circular ca se. The e lliptic 
Hill Hamiltonian is the followi ng JSzebehelv 19671; 
'ichtiaroglo ul Il98(]|: iMoons et a l. 1988; Llibrc & Pino] 
I199Q ; ,Brumberg fc Ivanoval Il990. ; Astakhov fc Farrellvi 
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12004 lAstakhov et~aLll2005l : iPalacian etal] 120061 ): 



H ^ E- 



(1+60 cos/©) 



3^ 



1 



8l3 



(1) 



Here E is the energy, r;, C^) = p defines the relative 
distance between the binary members mi and m2 and 
(PeiP'7iPc) = P is the corresponding momentum vector. 
The coordinate system rj, (") is rotating with constant 
angular velocity Jl© = (0, 0, 1) in the £, — rj plane. The 
eccentricity and true anomaly of the heliocentric orbit 
of the binary barycenter are e© and /q respectively. 
In this coordinate system the barycenter is located at 
the origin. The additive con stant is chosen such that 
the Lagrange saddle points l|Murrav fc Dermottlll999l ) 
in the circular (e© = 0) limit occur at i? = 0. 
Defining the reduced mass 1/2 

-=^^ (2) 

mi + m2 

allows the separate motions of mi and m2 to be recov- 
ered 



Pi 



-up 



P2 = (1 - i^)p- 

The equations of motion are the following; 



(3) 
(4) 



Ap3 



C 



(5) 



where A = 1/(1+60 cos /©) and dots denote derivatives 
with respect to time. In the circular Hill problem (6q = 
0) there exi sts an integral of the motion, the Jacob i 
constant, T l|Szebehelvlll967l : lMurrav fc Dermottlll999l '). 

r = 3e''-C'' + --(f + ?)'' + C'') + 8ii (6) 



2.2 Chaos in the Hill sphere 

If two bodies come within their mutual Hill sphere they 
cannot be captured into a binary permanently unless an 
energy loss mechanism is present. However, the partners 
can form a transient (quasi- bound) binary which may 
have a substantial lifetime (|Astakhov et all |2005| ). In 
order to form a quasi-bound binary the two primaries 
must first come close enough that mutual gravity dom- 
inates solar differential gravity. Thus, the proto-binary 



effectively defines a Hill sphere. At low energies the gate- 
ways to this region are the Lagrange saddle points, Li 
and L2, which separate interior and exterior (heliocen- 
tric) orbits. The interior region is the Hill sphere and 
has radius 



a© 



mi + m2 
3M, 







(7) 



Examination of Poi n care surfaces of section 
(iLichtenberg fc LiebermanI Il992|) in the Hil l prob 



lem |sim6 fc Stuchi 
lAstakhov fc Farrellvl 
or, equivalently. 




Astakhov et all 12003 : 
lAstakhov et al.1 l2005r) 

^ CRTBP for small masses 

(lAstakhov et al.ll2003h - reveals that, even at energies 
above the Lagrange saddle points, Li and L2, phase 
space is divided into regu lar Kolomogorov-Arnold- 
Mose r (KAM) regions (|Lichtenberg fc LiebermanI 
I1992I ). chaoti c regions and hy perbolic (direct scatter- 
ing) regions l|Zaslavskvl Il985h . Most importantly, the 
chaotic orbits separate the regular from the hyperbolic 
regions. Because incoming trajectories from outside 
the Hill sphere cannot penetrate the regular KAM 
tori at all in 2-dimensions (2D) and enter re gular 
regions exponentially slowly in 3D (|Nekhoros hcv 19771; 
Perrv fc Wigginslll994l : ILichtenberg fc LiebermanT992l : 
Astakhov et al.ll2003l '). particles entering the Hill sphere 
from heliocentric regions must either enter chaotic 
layers or scatter out of the Hill sphere promptly. Those 
that enter chaotic layers may remain temporarily 
trapped within the Hill sphere on time scales that are 
orders of magnitude larger than might be expected 
from estimates based on prompt transit of the Hill 
sphere. This dramatically increases their chances of 
being permanently captured by, e.g., gravitational 
scattering with a fourth body. 

It is interesting to note that chaos in the Hill sphere 
is similar to the situation in several problems of cur- 
rent interest in atomic and molecular physics; these in- 
clude the creation of non- disp ersive electronic Trojan 
wave p ackets in Rydberg atorn s (iBialvnicki-Birula et al.l 
19951 : iFarrellv fc Uze j Il995l : iLee et all Il995l . Il997l . 
200d : iKalinski et al.1 120051 ): the dynamics and ion- 



ization of Ry dberg atoms and molecules in exter- 
nal microwave llGriffiths fc Farrellvlll992l:lFarrellvlll994l: 



Farrellv fc Uzer"l995': 'Farrellv et al."l99g: 'Deprit et al.l 
1996: Brunello^et al. 1997; Bcllomo et al. 1997) or mag- 
netic fields fjohnson et al.lll983l : ISaini fc Farrellvll 19871 : 
lUzer et al.lll991 ); and the dynamics of ions in ion traps 



Howard fc Farrellvl [1993I : lElipe et al.ll2002l ). 



2.3 



Penetration of the mutual Hill sphere by 
the primaries 



The p roblem is similar to that studied bv lHenon fc PetitI 
(|l986l ) who investigated satellite e ncounters in the pla - 
nar circular Hill problem - see also lYoder et all l|l983h : 
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Impact parameter Eccentricity 




Figure 1. (colour online) Back-to-back histograms showing the number density of orbits entering the Hill sphere as a function 
of the parameters, b, e and i and the relative velocity v - see eq. I|12|l - in the circular [green (online), light grey (print) and 
elliptic [blue (online), dark grey (print)] Hill problems. These quantities are all defined in the asymptotic region. Approximately 
350,000 trajectories were integrated in both the circular and elliptic cases. The area under each histogram is normalized to 
unity; / is the probability density. All units are Hill units. The Jacobi constant is picked randomly as described in the text. 



IPetit &: HenonI l|l986l ): The two primaries are initially 
assumed to follow elliptical heliocentric orbits with 
semi-major axes a\ and 02 with some velocity disper- 
sion, V , around an elliptical Keplerian orbit lying in the 
invariant plane with semi-major axis Oq. Assume that 
the Keplerian velocity, Vk , is much greater than the rel- 
ative velocity of the two primaries as is thought to have 
been the case in the p rimordial Kuiper-belt, i.e., prior to 
dynamical excitation ((Chiang et"aL 20061: Levison et al.l 



l2006l : lGladman fc ChanlbOOd : iQuillen et al.l l2004'): How 

do the values of their orbital elements "at infinity" de- 
termine how close the primaries will approach and, if 
they do so approach, then is the encounter favourable 
for binary production? To answer these questions we 
first need to consider the asymptotic pre-encounter be- 
haviour of the system. 



action potential. Integrations originating in the asymp- 
totic regime must start at some finite t = to < 00. 
In principle, one could start the integrations at ex- 
tremely large separations such that the mutual inter- 
action terms in eqs. (O are sufficiently small. Alter- 
natively, one can use asymptotic expansions to start 
the incoming solution (from infinity) and to continue 
the outgoing solution (to infinity). In between the so- 
lution is obtained numerically by integrating the equa- 
tions of motion. High-order asymptotic solutions have 
been derived by Henon and Petit in the planar circu- 
lar H iU prob lem Jllenon fc PetitI Il986l: iPetit fc HenonI 
1 19861 ) and bv lBrumberg fc Ivanoval l 1990l ) in the three- 
dimensional elliptic Hill problem. The asymptotic solu- 



2.3.1 HiU orbital elements 

The asymptotic behaviour, when the relative separation 
between the primaries, p, is large {\t\ — » 00), is compli- 
cated by the infinite range nature of the "1/p" inter- 
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tions are, to lowest order, 



e = 



s{t - 



rj = ~—a{t — <j)) + 2esm{t — T) 

C, — i sm{t — Lu) 
^ = e sin(t — r) 
3 

■n = 

= i cos(t — u) 



- a + 2 e cos(t — r) 



(8) 



where a, e, and i are sometimes called the Hill, or "re- 
duced," orbital elements: here b — \a\ is the impact pa- 
rameter, e is the eccentricity and i is the inclination in 
Hill units while r, 0, and uj are phase angles. Alterna- 
tively, b can be thought of as the fractional distance of 
the semi-major axis from the reference orbit a©. Ex- 
plicitly the reduced elements are related to the usual 
semi-major axis (flc), eccentric ity (se) and inclination 
(ic) i n the CRTBP as follows (jGreenzweig fc Lissauerl 
I1990I I. 

K - 1) 

a = ^ — 

Rh 

— 

" Rh 

^ = ^ (9) 

The set of Hill orbital parameters has been 
used extensively in studies of the accretion of plan- 
etesimals or the dynamics of pa rt icles in plane- 
tary rings ijNakazawa fc Idal 19881: Nakazawa et al] 



1989"; 'Ohtsuki fc Ida' '1990^; 'Wetherill fc Stewart' 'l989l; 
Grccnzwcig & Lissauor 1990, 1992; Wetherill fc St ewart 
19931: lOhtsuki et alj I2OO2I: Istewart fc Idal l200o[ 



Petit fc Henonlll987al lbt lRafikovll200ll . I2003allbll3 1dh 



In the circular limit the quantity (j> can be elim- 
inated by an appropriate choice of the origin of time 
resulting in five elements. Equations (|8} are for the rel- 
ative motion of the primaries; by attaching subscripts, 
j = 1,2, one may obtain corresponding expressions for 
the individual primary orbits. 

The Jacobi constant can be expressed in closed 
form in terms of the set a,e,i 



(10) 



For our purposes it is sufficient to use the lowest-order 
versions of the asymptotic expansions - eq. (|8]) - which 
correspond to Keplerian motion and which are valid in 
both the pre- and post-encounter asymptotic regions. 
Strictly, higher order corrections should be included; 
however, even if the interaction is not negligible the or- 
bit can, instantaneously, still be described by elements 
having the form of eq. ([8]). Assume that, in this way, 
we compute a set of orbital elements {a',e',i') which 



are approximations to the true orbital elements (a, e, i). 
Provided that the error is small then one can think of 
{a',e',i') as being the exact orbital elements but for 
a slightly different value of the Jacobi constant, F'. In 
our Monte Carlo simulations F is chosen randomly and, 
therefore, neglecting to use the higher order asymptotic 
formulae is not expected to affect the results. The main 
reason for not using higher-order formulae is the addi- 
tional computational effort involved in solving for the 
orbital elements, which must be done iteratively. 



2.3.2 The Hill velocity and the relative velocity 

Physically, the Hill velocity is the orbital velocity 
around a large body at the Hill sphere ra dius assuming 
no solar perturbations (iGoldreich et al.l ,2004; ,Rafikovl 
l2003dl : lMurrav-Clav fc Chiang||2006l ) - so that one revo- 
lution at the Hill radius completes in one orbital period 
around the Sun. However, this turns out to be a fairly 
approximate quantity because the dynamics is essen- 
tially three-body in the case of TNB formation. We will 
modify this definition slightly to include the case where 
the Hill sphere is defined by a binary, rather than by a 
single mass, but for which the separation between part- 
ners is considerably less than Rh, 



Vh 



G(mi + 7712) 



Rh 



Rh. 



(11) 



In HiU units, the Hill velocity vh = 3"^'''^ ~ 1.2. To see 
the correspondence in physical unit s, consider a binary 
with the following characteristics (jVeillet et al ]|2002'): 
(a) a barycenter semi-major axis of 45 AU; (b) 100 km 
radii binary partners; and, (c) density d = 1 g cm~^. 
For these parameters Rh ~ 7.5 x 10^ km, mi = 1112 ~ 
4.0 X lO^^g, the Keplerian velocity, Vk ~ 4.4 km/s and 
the Hill velocity Vh ~ 0.9 m/s. 

Scattering is said to be in the shear-dominated 
regime (dynamically cold) when v < vh and in the 
dispersion-dominated regime (dynamically hot) if u > 
VH (S tewart fc Id a 2000; Collins fc Sari 2006). Shear 
is induced by the difference between the Keplerian 
angular velocities of primaries having different helio- 
centric semi-major axes. If the relative approach ve- 
locity of two particles is greater than the differential 
shear in the disk across the Hill (or ti dal) radius then 
the dynamics is disp ersion-dominated (jColdreich et al.l 
I2OO4I : iRafikwll2003cl ) The parameter which determines 
whether the system is hot or cold is the ratio of the ve- 
locity dispersion to the shear across the Hill radius, i.e., 
V/Q.Q Rh = v/vH iRafikov 2001). Thus, the HiU veloc- 
ity serves roughly to demarcate the transition from two- 
to three-body dynamics for encounters between "small" 
{R < 1 km) and "big" {R > 100 km) bodies; i.e., ff 
the relative speed of bodies undergoing close encoun- 
ters is greater than the Hill velocity then two-body dy- 
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Circular Hill Problem Elliptical Hill Problem 
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Figure 2. (colour online) Kernel smoothed scatterplots showing log]^Q(THiii) where Tniii is the Hill lifetime (in years) in the 
circular (left) and elliptic (right) cases as a function of the relative asymptotic velocity v scaled by vjj - see eq. II12I I. In Hill 
units 300 years ~ 27r; i.e., the Hill unit of time ~ 50 years. Scale runs from low [yellow (online), light grey (print)] to high 
[red (online), dark grey (print)]. In this representation smoothed point densities are computed using a k ernel density estimator . 
In regions of sparse data the plot is augmented by black dots which represent individual data points l lGentleman et aI1l200d 
ICarr et al.lll98'7lMGentlemanll2007^ . For further discussion of kernel smoothing see Appendix A. 



namics is expected to provide a good approximation to 
the dynamics otherwise solar tides mu st b e included, 
i.e., three-body effec ts are important (jRafikovi i2003dl : 
iGoldreich et al.ll2004 ). Hovy ever, this is not a sharp tran- 
sition and iRafikovl (|2003dh has identified a transition 
regime. In this work we find that binary capture spans 
such a transition region. 

iRafikovl (|2003dh has further argued that the growth 
of big bodies (planetary embryos) in a planetesimal disk 
undergoes a direct transition from a relatively short run- 
away growth phase to a much longer oligarchic phase. If 
TNBs formed during the longer oligarchic phase then 
this would suggest that the velocity dispersion was 
larger than vh- In contrast. IGoldreich et al.Tl|2002l ) have 
estimated that the velocities of ~ 100 km sized bod- 
ies were on the ord er of fg/S in t he ea rly Kuiper-belt. 
As pointed out by IChiang et al.l (|2006l ) if the velocity 
dispersion, v, of big bodies is less than vh then these 
bodies will collapse into an essentially two-dimensional 
disk due to so-called runaway cooling induced by the 
dynamical friction exerted by the sea of small bodies. 
Therefore, TNB mutual orbit inclinations ought to be 
similar, i.e., mutual orbit normals will be approximately 
parallel. In fact, TNB inclinations appear to b e ran- 
domly distributed (|No11||2003| : IChiang et al.ll2006l ) . This 
implies an isotropic primordial velocity distribution and 
suggests that t he big bodies did not all originally lie in 
the same plane (|No11 et al.ll2007h and that TNBs formed 
during the oligarchic phase. 

While these arguments lead to the inference that, in 
the early Kuiper-belt, v > vh and that observed TNB 
inclinations are primordial, it could also be the case 
that the observed TNB inclinations are the result of ei- 
ther (i) dynamical stirring of the Kuiper-belt after bina- 



ries formed jChiang et al.|[2006l : iLevison et al. 20061) or 
(ii) post-capture bin ary hardening (jHeggie et al T fl996l : 
iHeggie fc Hutll2003l ). or both. In our previous simula- 
tions we found that inclinations do no t change signif- 
icant ly during the hardening process (lAstakhov et al.l 
l2005l ) - this is similar to the approximate conserva- 
tion of irregular sate llite inclinations unde rgoing CAC 
at Jupiter or Saturn (|Astakhov et al.ll2003l ). 

A related consideration is deciding how the "rel- 
ative velocity" is best defined. As the two bodies ap- 
proach each other from infinity then their relative ve- 
locity changes and, if the bodies enter a chaotic layer 
within the Hill sphere, fluctuations in the relative ve- 
locity can be large. At large separations we adopt the 
following definition, 



Rh 

a© 



(12) 



This expression is similar to the average velocity of a 
planetesimal, relative to other planetesimals in a swarm 
with mean Keplerian velocity vk, and as obtained by 
averaging over an epicycle and a vertical oscillation 
(jCreenzweig fc LissaueJll990l . Il992l ). It turns out that 
the relative velocities at infinity which lead to penetra- 
tion of the Hill sphere - see Fig. [1]- are only roughly on 
the order of the Hill velocity. 



2.4 Simulations of primary encounters 

The most straightforward approach might appear to be 
a Monte Carlo simulation in which the equations of rel- 
ative motion - eq. ([5} - are integrated inwards from large 
initial separations. However, this is not practical com- 
putationally because most initial conditions picked ran- 
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domly at infinity will not lead to primary encounters, 
i.e., separations less than Rh- Especially is this true 
for the three-dimensional elliptical problem. Instead we 
adopt a different approach which involves the several 
stages now to be described. 

2.4-1 The feeding zone 

The first task is to try to establish what ranges of initial 
conditions at infinity can result in penetration of the Hill 
sphere. If these ranges can be pinned down they define 
a feeding zone in phase space. Assume that the edges of 
the presumed zone can be approximately delimited by 
maximum and minimum values of the three orbital ele- 
ments. The most obvious constraint is that the Jacobi 
constant should have a value higher than its value at 
the Lagrange points, i.e., F > 0. We desire to obtain, in 
addition, constraints on the ranges of each of the orbital 
elements if such constraints exist. Of course, the mere 
fact that a set of orbital elements is contained within 
the feeding zone need not imply that the corresponding 
trajectory will actually enter the Hill sphere. It is also 
possible that initial conditions lying outside such a feed- 
ing zone will enter the Hill sphere. The feeding zone is, 
therefore, an approximate, but potentially useful, con- 
cept because it narrows down the ranges of the orbital 
elements which can, in principle, lead to long-lived en- 
counters within the Hill sphere. Reasonably good con- 
straints on the feeding zone can be found numerically 
as follows. 

(i) Generate initial conditions randomly inside the 
Hill sphere for values of F chosen randomly and uni- 
formly in the range F £ (0, 8) in Hill units. 

(ii) Cartesian positions and velocities are then gener- 
ated uniformly and randomly inside a sphere of radius 
1.2Rh. 

(iii) Next Hill's equations are integrated until the tra- 
jectory penetrates the surface of a sphere, ^2, in the 
asymptotic region; here S2 is chosen to have radius 
~ 3507?£f. At this point the integration is stopped and 
the orbital elements at infinity stored. 

(iv) Some initial conditions will, of course, have 
been started within (and, in three-dimensions, between 
(|Nekhoroshevll 19771 )) regular KAM zones inside the Hill 
sphere. These initial conditions are of no interest for 
present purposes because they lie in regions that cannot 
be penetrated at all in two-dimensions and only expo- 
nentially slowly in higher dimensions. Therefore, they 
are discarded; however, first they must be identified. 
This is accomplished as follows: if, after a sufficiently 
long time, Tcut, a trajectory has not passed through 
S2 then that orbit is discarded. Of course, it is possi- 
ble that discarded trajectories did not actually start in- 
side KAM regions but, instead, were trapped within ex- 
tremely long-lived chaotic zones - i.e., they are amongst 



the initial conditions of greatest interest. Therefore it 
is important that T^ut be chosen large enough that (a) 
such cases are relatively few and (b) the results are in- 
sensitive to its precise value. We chose T^ut = 1000 in 
dimensionless Hill units or ~ 50, 000 years. While this is 
shorter than the very longest lifetimes shown in Fig. [2] 
this method can still detect such orbits because it is un- 
likely, in this approach, that initial conditions selected 
at random inside the Hill sphere will subsequently spend 
their entire Hill lifetime inside this region. 

(v) Because we integrate a very large number of tra- 
jectories we are confident that this approach allows for 
the harvesting of essentially all types of initial condi- 
tion which, coming in from infinity, will penetrate the 
Hill sphere. 

Figure [T] compares histograms of the asymptotic 
distributions of orbital elements and the velocity of ini- 
tial conditions, which, started in the asymptotic regime, 
go on to penetrate the Hill sphere in the circular and 
elliptical problems. The distributions in the circular and 
elliptical cases are generally similar except that the el- 
liptical distributions are fatter; this is most marked in 
the distribution of mutual orbit eccentricity and veloc- 
ity. In these simulations the heliocentric eccentricity was 
confined to the range £ (0,0.3). 

Roughly speaking, in the circular Hill problem, the 
feeding zone is defined by the following ranges (in Hill 
units): r e (0,4),& e (0.8, 5), e € (0,5),|i| € (0,1.5) 
and V £ (0, 4). These ranges are approximate and initial 
conditions lying outside them may lead to trajectories 
which penetrate the Hill sphere but, we will find, they 
tend to do so promptly, that is they do not get caught 
up in chaotic layers. Similarly, not all initial conditions 
lying inside these ranges necessarily pass through - or 
close to - the Hill sphere. 

2. 4-. 2 Quasi-binary Hill lifetimes 

We define the Hill lifetime, Tm\i, for each set of initial 
conditions at infinity to be the total time the resulting 
trajectory spends within the Hill sphere - see Table 1. 
To calculate Tniii the trajectories used in identifying the 
feeding-zone are back-integrated from infinity and their 
time inside the Hill sphere recorded. 

The key to the CAC scenario is the dramatic ex- 
tension in the time the binary partners spend within 
the Hill sphere due to their having become entangled 
in very long lived - though chaotic - regions of phase 
space. We computed lifetimes of all trajectories which, 
starting in the asymptotic region, end up penetrating 
the Hill sphere. The results are shown in Fig. [2] where 
Hill lifetimes are plotted as a function of the relative 
velocity at infinity obtained using eq. (|12p . The large 
island lying between « 10 and « 100 years corresponds 
to essentially direct or "ballistic" transit through the 
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Hill sphere and appears in both the circular and elliptic 
cases. In both cases a second large island exists and cor- 
responds to lifetimes on the order of thousands of years 
or more with much longer lifetimes also being possible. 

Figure [5] can be thought of as a nonlinear map 
which relates asymptotic velocity to transit time 
through the Hill sphere. Clearly, two (or more) tra- 
jectories with the same asymptotic velocity can spend 
very different times inside the Hill sphere. This map also 
shows that very small asymptotic velocities tend either 
not to enter the Hill sphere at all, or, if they do enter, 
they transit rather quickly. In part, this is because the 
actual velocity after entering the Hill sphere can be con- 
siderably greater than the asymptotic value. Note that 
Fig. [2] covers essentially the full spectrum of asymptotic 
velocities which can lead to Hill sphere penetration. 

The elliptic case shown in Fig. [2Ib) involves a much 
larger range of asymptotic velocities and orbital ele- 
ments than does the circular problem. Simulating the 
four-body dynamics - i.e., intruder - binary scattering 
inside the Hill sphere - in the elliptical case will be sig- 
nificantly more demanding computationally than in the 
circular case because of the larger space of intruder and 
binary initial conditions, e.g., both the primaries and 
the intruder are now allowed to follow elliptical orbits. 
Mainly for this reason, from now on, in this paper, we 
specialize to the circular Hill problem but consider the 
elliptic problem worthy of further study. 



3 STABILIZATION BY INTRUDER 

SCATTERING IN THE CIRCULAR HILL 
PROBLEM 

Having formed a transient binary the next step is to 
capture it permanently. In this section we investigate 
how the masses and velocities of the primaries and the 
intruder affect nascent (i.e., just-captured) binary for- 
mation as well as the properties of any binaries that 
result. The simulations are done in the four-body Hill 
approximation (,Schccrcs.l998, ') as is now described. 



3.1 Four-body Hill approximation 

Three comparatively small bodies, with a mutual 
centre-of-mass, Rc, orbit a much larger body - the Sun 
- mo = 1 on a near-circular orbit e.g., the primaries 
and the intruder. The total mass of the three bodies is 
defined by 



^m, « 1 



(13) 



where Rc ~ a = (1,0,0) defines the motion of the 
three-body centre-of-mass along an almost circular orbit 



which defines the rotating fra me. The vector equations 
of motion are (|Scheereslll9"9^ ) 



p + X [2p + X p] = -p + 3a(a ■ p) - (ai + 02)7^ 



+ as 



Pi - P2 P3 - Pi 



(14) 



Pi + Vlx [2p3 -I- n X ps] = -p3 + 3a(a ■ ps) 



P3 - Pi P3 - P2 



027 



P3-piP IP3-P2P 



(15) 



Here ps is the coordinate of the third intruder body, 
m-i, and nij = jj,aj where 



(16) 



When ma — eg. II 1411 reduces to the 
three-body Hill problein dffill Il878l : ISzebehelvl 1 19671 : 
iMurrav fc Dermotd[l999h and becomes uncoupled from 
eq. (fT5)) . 



3.2 Numerical procedure 

This numerical simulations were performed as follows. 

(i) Initial conditions for the primaries "at infinity," 
and which are guaranteed to penetrate their mutual 
Hill sphere, are generated as described in Subsec. 2.4. 
That is, initial conditions are generated randomly inside 
the Hill sphere and integrated until the trajectory pen- 
etrated a sphere, S2 , of radius ~ 350Rh at which point 
the integration is stopped and the orbital elements com- 
puted. 

(ii) The integration is then run backwards so as to 
compute the Hill lifetime of the orbit and the orbital 
elements are stored. 

(iii) This procedure is repeated until a sufficiently 
large cohort of initial conditions at infinity has been 
generated. 

(iv) Intruder initial conditions at infinity are gen- 
erated in a similar way except that (a) initial con- 
ditions for the outward integration are generated in- 
side a sphere of radius 2Rh and (b) the radius of the 
sphere "at infinity" was chosen randomly in the range 
30Rh < S2 < 350Rh. This is done to allow for the pos- 
sibility of stabilization by "near-misses", i.e., intruders 
which do not actually penetrate Rh and also to ensure 
that the phases of the intruders with respect to the bi- 
nary orbit are varied. 
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Figure 3. Histograms of Hill lifetimes in the pilot calculations described in the text. Frame (a) shows the relative frequency of 
orbits with Hill lifetimes in the range Thhi £ (0, 1000) years. Frame (b) shows the frequency of captured orbits as a function of 
the Hill lifetime. In computing the histogram in (b) equal numbers of initial conditions were used in each 100-year interval as 
described in the text. 



(v) Next, the primaries are integrated back from in- 
finity in the three-body HiU approximation until they 
come within a radius Rh < Ro < 5Rh of each other. 
The actual radius, Ro, is generated randomly and the 
orbital parameters of the primaries are then stored. 

(vi) The mass ratio of the primaries is selected ran- 
domly and the stored orbital parameters of the pri- 
maries are rescaled using the actual values selected for 
mi and m^. This is possible because, in the three-body 
Hill approxi mation, all masses can be scaled out of 
the problem l|Murrav fc DermottI Il999l ). However, the 
masses appear explicitly in the four-body Hill problem 

- see eq. (fH)) . 

(vii) The mass of the intruder is selected randomly up 
to the mass of the larger primary. Thus the maximum 
ratio of intruder mass to total binary mass is 0.5. 

(viii) The intruder is then launched from infinity and 
the integration of the primaries simultaneously resumed 

- this time using the four-body equations of motion. The 
system is integrated until either the binary breaks up or 
it is stabilized. We have confirmed numerically that at 
the start of the four-body integrations the primaries and 
the intruder are sufficiently well separated that using 
re-scaled three-body Hill initial conditions to start the 
four-body integrations is legitimate. Thus, these compu- 
tations cannot describe simultaneous "three-body plus 
Sun" encounters in which the primaries and the intruder 
all interpenetrate their mutual Hill sphere at about the 
same time. Of course, these encounters are describable 
by the four-body Hill equations but they are extremely 
rare and so are neglected. 

(ix) Binary stabilization is registered if the binary 



survives for 10 times longer than its lifetime inside the 
Hill sphere in the absence of intruder scattering, or 
for T = 200 Hill units (~ 9640 years), whichever is 
larger. For example, if a transient had a Hill lifetime 
greater than 200 Hill units, say 201 Hill units or ~ 
9,693 years years then we counted it as stabilized only 
if it survived for at least 96,930 years. Capture statis- 
tics were quite insensitive to using lifetime extension 
multipliers larger than 10 but were somewhat sensitive 
to using multipliers < 7. As expected, intruders can 
prolong quasi-binary lifetimes by pushing trajectories 
deeper into chaotic zones but without actually stabiliz- 
ing them permanently. 

(x) Keeping all masses fixed this procedure is then 
repeated for 5 different values of Rq. This has the effect 
of sending the intruder towards the binary at different 
relative configurations of the binary partners. 

(xi) For each binary orbit the overall process was re- 
peated for 1000 randomly selected intruders, each time 
varying all masses randomly. 

(xii) Individual integrations were stopped and results 
discarded if particles came within a distance ta = 10~^ 
Hill units of each other. This radius is somewhat arbi- 
trary but is roughly the radius of a typically-sized binary 
partner and thus this choice corresponds approximately 
to a physical collision between the binary members. Col- 
lisional singularities could be a voided by regularisation 
(|Szebehelvil"967l : I Aarsethll2003l ) but we preferred to stop 
the integrations if a collision occurred. 

(xiii) In total 15,000 quasi-bound binary orbits were 
harvested. For each binary, 1000 randomly selected in- 
truders from a pool of 15,000 were then sent towards 
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it as described above. Each binary-intruder encounter 
that led to capture was considered to be a single cap- 
ture event. 



3.2.1 Pilot calculations to identify a threshold Hill 
lifetime 

Examination of Fig. [IJa) reveals a potential computa- 
tional difficulty associated with the algorithm just de- 
scribed. Ballistic trajectories - i.e., those which have Hill 
lifetimes ~ 50 — f 00 years and follow hyperbolic or near- 
hyperbolic orbits - significantly outnumber trajectories 
which penetrate and become entangled in chaotic lay- 
ers. Because, as we will show, these trajectories have 
very low capture probabilities, including them in the 
full simulations would swamp the calculations. On the 
other hand, we must first demonstrate that such capture 
probabilities are low. This is als o important because the 
TNB formation mechanism of iGoldreich et al.l (|2002l ) 
assumes that it is precisely these trajectories which lead 
to capture, i.e., trajectories which transit the Hill sphere 
on time scales on the order of Rh/Vh = l/f^Q ~ 50 
years. 

We therefore performed a set of pilot calculations 
in which the capture of binaries with Hill lifetimes in 
a limited interval, i.e, Tuiii £ (0, 1000), years was stud- 
ied. Figure |3ja) shows the initial distribution of Hill 
lifetimes obtained in this range using the procedure de- 
scribed above. It is apparent that if this distribution 
of orbits were to be used directly then the vast major- 
ity of the computations would involve trajectories with 
Tniii < 100 years. In fact, simulations using this distri- 
bution had to be abandoned because most of the inte- 
grations of these initial conditions did not lead to cap- 
ture. However, even in these limited simulations it was 
clear that the capture probabilities of trajectories with 
longer Hill lifetimes dominated despite their relative un- 
der representation in Fig. [3ja) . 

To study capture probabilities further we developed 
an alternative procedure. 

(i) 100 initial conditions were selected randomly in 
the ten intervals: Tniii G (0, 100) years, Tmn G (100, 200) 
years, etc., up to THiii G (900, 1000). 

(ii) 2000 intruders - selected as described above - were 
sent in towards these proto-binaries and the number of 
scattering events which led to capture in each interval 
was recorded. 

The results are shown in Fig. [Sfb). The probabil- 
ity of capture for THiii < 100 years is extremely small. 
Of course, in the original disk, primaries which enter 
their mutual Hill sphere and escape without being cap- 
tured can, at a later time, re-enter the Hill sphere and 
new opportunity to become caught-up in a chaotic layer 
presents itself. That is, escape is from the mutual Hill 



sphere and not from the disk itself. Therefore, even 
though the probability of a ballistic binary being cap- 
tured in a single pass might be small repeated passes 
through the Hill sphere are possible. The overall cap- 
ture probability will clearly depend on the mass and 
orbital element distributions in the original disk. This 
is not included in our simulations. 

Based on Fig. |3] we included only binaries with 
Tmu ^ 200 years in the large set of simulations now 
to be described. We note that the simulations in the 
pilot calculations and the full calculations show similar 
dependencies on intruder mass and binary mass ratio. 



4 RESULTS 

4.1 Effect of intruder mass 

The masses of the primaries were varied randomly and 
the mass of the intruder was then chosen randomly up 
to the mass of the larger of the two primaries. Figure 
Ua) indicates that the efficiency of binary capture falls 
off with increasing intruder mass, i.e., intruders of com- 
parable size to the primaries tend to (i) leave the proto- 
binary essentially unaffected, (ii) destabilize the rather 
delicate quasi-bound binary (i.e., reduce its natural Hill 
lifetime by causing the complex to break up before its 
time) or (iii) stabilize it against ionization with lower 
probability than do smaller intruders. Destabilization 
happens, e.g., when a three-body resonance forms which 
has the effect of causing the binary to split up prema- 
turely. However, we have also found some rarer cases 
where the intruder forms a resonance which actual lives 
longer than the Hill lifetime but which, nevertheless, 
eventually breaks up. 

The rate at which the tail of the distribution (large 
intruders) falls off depends to some degree on the 
closest-approach distance of the intruder to the binary 
- see Fig. (5) That is, large intruders can be effective sta- 
bilizers provided that they do not get too close to the 
binary. On average, small intruders can penetrate much 
deeper and still stabilize the binary. It is apparent from 
Fig. |4ja) that small intruders are most efficient at sta- 
bilization. Given that roughly equal diameter binaries 
dominate in Fig. Iljb) we conclude that small intrud- 
ers tend to lead to roughly equal mass ratios. There- 
fore, in view of Fig. [S] small Rrain values also correlate 
with roughly equal mass ratios. The reason for this ef- 
fect is that small intruders can penetrate deeper before 
causing the binary to break up. For unequal mass bi- 
nary partners this means that they can more effectively 
get caught-up in three-body resonances which tend to 
destabilize the binary. 
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Figure 4. Histograms showing the probabihty density - i.e., the area under each histogram is normahzed to unity - of per- 
manently captured binaries as a function of: (a) intruder mass where rrtp = m3/(mi + 012) X 100 is the mass of the intruder 
expressed as a percentage of the total binary mass; (b) diameter ratio, D2/D1, of the primaries; (c) computed visual magnitude 
difference, Amag', (d) logj^QAmag; (e) average angular momentum, /lavo = prograde orbits correspond to /i^ > and 

retrograde orbits to < 0, and; (f) logj^g of their Hill lifetime in years, i.e., the maximum time the quasi-binary would spend 
inside the Hill sphere in the absence of intruder scattering. 
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Figure 5. Scatter plot showing the minimum penetration distance, Rmin, of intruders which successfully capture binaries, as a 
function of intruder mass expressed as a percentage of the total binary mass. The attendant histograms show the distributions 
of Rmin and intruder mass. 



4.2 Effect of binary mass ratio 

As noted, the binary mass ratio was also varied ran- 
domly throughout the simulations. Assuming that the 
binary partners have the same density, mass ratios can 
then be translated into size- or diameter-ratios. Figure 
Ub) shows a strong preference for roughly same-sized 
binary partners although this might be offset if a more 
realistic, e.g., a power-law, mass distribution were used. 



4.3 Magnitude differences 

Observed magnitude differences, Amag, of binaries 
can be used to obtain informati on about the rela- 
tive sizes of the partners involved jTruiiUo et al.ll200ll : 
INoU et all l2006al : iPetit et all 120061 ). Assuming that 
the primari es have the same density and albedo as 
each other ('Cruikshank' '2005'; 'Lvkawa & Mukai| l2005l : 
[Stansbcrry ot a l. 2006; Petit & Mousis 2004) allows a 
relationship to be established between binary mass ra- 
tios obtained from calculations and observed Amag val- 
ues. The diameters o f the primaries, Di, D2, are related 
to A mag as foUows (jsheppard fc JewittI |2002| : iHughesI 
I2OO3I '): 



where, by (our) definition, D2 < Di . Figures EJc) and 
(d) show the predicted distributions of Amag, and (so as 
to amplify the region around Amag = 0) logj^f^(Amag)- 

In order to ma ke a more direct connection with 
actual observations (|No11 et al.l l2006al 'l it is important 
to know not only the magnitude differences but the 
magnitudes themselves. Our model makes no predic- 
tions about the absolute sizes of the binary partners - 
only their mass ratio is predicted. However, various fits 
have been made to the size distribution of TNOs with 
some suggestio n that the number of objects with radii 

< 40 - 70 km (iBernstein et al.ll2004 I2OO6I : |Pan fc Saril 
J' „_ I i^T-- '11 11 ' . " , \ fj L . : — 



£1 



= 10" 



(17) 



I2OO5I : lEUiot et all I2OO5': Petit ct al. 2006') is somewhat 
less than expected. LBcrristoin ot al. ( 2 004i ) have pro- 
posed a double power-law fit to the differe ntial surface 
densit y of TNOs as a function of magnitude. iPetit et al.l 
(|2006l ) have argued that a double power-law expan- 
sion is unnecessary although they provide parameters 
which fit their observa tions to the double power-law of 
IBernstein et all (|2004| ). For these reasons we used the 
single power-law distribution of iPetit et al.l l|2006h to 
simulate a plot of magnitude against magnitude differ- 
ences. This was done as follows; 

(i) The differential magnitude distributio n was mul- 
tiplied by the efficiency function reported bv lPetit et al.l 
(|2006l ). which is a product of two hyper boUc tangents. 
This produced a "corre cted" distribution (|TruuUo et al.l 
I2OOII : lEUiot et al.ll2005h . 
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(ii) This distribution was normalized according to the 
number of captured orbits obtained from the simula- 
tions (~ 70,000). 

(iii) Using this distribution visual magnitudes were 
randomly assigned to the larger member of each nascent 
binary with magnitudes, m, in the range 20 ^ m ^ 30. 

(iv) A bootstrapping procedure was used in which 
the previous step was repeated 5 times using a different 
random number seed in each case. Thus, each captured 
binary actually enters in with 5 different randomly as- 
signed magnitudes. This procedure improves the statis- 
tics for smaller magnitudes where comparatively few in- 
stances are generated. 

The results are shown in Fig. |6] Keeping in mind 
the caveats already noted in regard to albedo and den- 
sity, the results in Fig.|3]and Fig.[S]demonstrate that the 
capture model predic ts a similar clusteri n g of A .maa < 1 
to that observed by INoU et al] (|2006bl . l2007h . In our 
model, the physical reason is chaotic destabilization of 
unequal mass quasi-binari es in comparison to tho se hav- 
ing roughly equal masses l|Astakhov et al.ll2005^ . 

4.4 Binary orbit mutual inclinations 

Strictly speaking it is not possible to compute meaning- 
ful orbital elements for nascent, i.e., just-formed, bina- 
ries. This is because th ese objects gener ally follow highly 
non-Keplerian orbits l|Astakhov et al.l [2005) and their 
orbital elements will change due to post-formation hard- 
ening. However, in previous studies of the CAC model 
we have found that the sign of the angular momen- 
tum component, hc_ = ^Pn — rjpe^, tends to be roughly 
conserved and that the inclination is quite robust to 
post-capture hardening and orbit reduction though gas- 
drag jA stakhov et al. 200j|) or gravitational scattering 
(jAstakhov et al.ll2005i '). However, the inclination distri- 
bution may be sensitive to the hardening mechanism as 
well as other assumptions in the model. Further inves- 
tigation is called for. 

During the transition from proto-binary (slightly 
unbound) to nascent binary (slightly bound) the incli- 
nation hardly varies during the rather "delicate" inter- 
action with an intruder (this is not always the case, 
but it is the case more often than not). Thus the incli- 
nation distributions before and after capture are quite 
similar. This is a direct consequence and signature of 
CAC because it is the location of the chaotic layers 
in phase space that determines which inclinations are 
most favourabl e for a given value of the J acobi con- 
stant (Astakhoy et a,l.ll2003l : lAstakhov fc Farrellv 2004 : 
lAstakhov et al.ll2005l ). 

For each nascent binary we computed its average 
post-capture angular momentum, /lavc = (^c)- These 
results are shown in Fig. |4je). The distributions are 
quite asymmetric, especially in the tails, and predict a 



broad distribution of orbit inclinations. Although diffi- 
cult to determine observationally, this suggest that mu- 
tual inclination distributions will provide a strong test 
of our model. As noted, this result is similar to the CAC 
model of irregular satellite capture which predicts differ- 
ent sub-populations of retrograde and prograde moons 
at Ju piter and Saturn (|Astakhov et al.ll2003l : iHamiltonl 
|2003| V However, further study of the details of the ef- 
fect of hardening on inclination distributions is clearly 
merited. 



5 COMPARISON WITH THE L^S AND 
MODELS 

iGoldreich et~aLl (|2002l ) estimated binary formation rates 
corresponding to two different mechanisms. The first 
mechanism, denoted "L^s," involves dynamical friction 
where "L" refers to a "big" body and "s" to a sea of 
small bodies which contain most of the mass. Essen- 
tially, two big bodies (L -|- L) enter their mutual Hill 
sphere and become permanently captured by the dy- 
namical friction exerted by the sea of small bodies. The 
second mechanism, denoted "L^," involves gravitational 
scattering by a third big body, i.e., 

L^*+L^L^+L* (18) 

which extracts energy (*) from the proto-binary thereby 
stabilizing it. 

5.1 The L'^s mechanism 

Dynamical friction, in principle, provides an efficient 
route to binary hardening and orbit reduction. How- 
ever, it is not sensitive to the mass ratio of the initial 
quasi-bound binary. The Hill equations, even including 
dissipation, can b e scaled so as to contain no par ameters 
(|Astakhov et al.ii2003. : .Murrav fc Dermottlll999l ) and so 
do not depend on the masses of the binary partners. 
Thus, no mass selection is predicted. In addition, esti- 
mates of planetesimal mass distributions suggest that 
L^s m ight not be an efficient binary production mecha- 
nism (jFunato et all [20041 '). 

5.2 The LP' mechanism 

In L"^, as originally proposed, roughly equal mass bina- 
ries result directly from the assumption that a transient 
binary undergoes a scattering encounter with a third 
body of comparable mass to the binary partners (which 
are also assumed to be of similar mass) . In essence, the 
preference for order unity mass ratios is partly a result 
of the assumptions of the model. Furthermore, multiple 
encounters are assumed to harden the binary. We 
have shown previously that multi ple large body encoun - 
ters tend to disrupt the binary (|Astakhov et al.l^2005^ . 
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Figure 6. (colour online) Kernel smoothed colour density scatterplot showing the visual magnitude (mag.) of the larger binary 
primary as a function of magnitude differ ence, Amag- T he distribution used to assign visual magnitudes is modified by a double 
hyperbolic tangent efficiency function (.Petit et al.ll2006h . A total of « 70,000 captured binary orbits were used over a magnitude 
range 20 ^ m ^ 30. Colour scale runs from low [yellow (online), light grey (print)] to high [dark blue (online), dark grey (print)] 
and the kernel smoothing is as described in Fig. |2] and Appendix A. 



Examination of the actual orbital elements produced in 
these simulations indicate that neither a single nor a 
small number of L"^-type events is sufficient to produce 
Keplerian orbits. Those simulations - together with the 
current set - suggest that the assumption that intruder 
scattering stabilizes binaries with a probability of order 
unity is a considerable overestimate, especially for bina- 
ries which pass through the Hill sphere on time scales 
of order Rh/Vh- 

However, both o f the mechanisms presented by 
iGoldreich et al.l (|2002l ) yield binary formation rates on 
the order of 10~® or 10"'^ per year which are compara- 
ble to large body growth time scales and are, therefore, 
quite reasonable in themselves. However, these numbers 
depend sensitively on assumptions made about various 
quantities, e.g., the residence time of a transient binary 
inside the Hill sphere and, in the case of L^, the prob- 
ability of a fourth body entering the Hill sphere and 
stabilizing the transient binary during its Hill lifetime. 
Getting numerically reasonable rates does not necessar- 
ily imply that the underlying mechanism and assump- 
tions are reasonable. In fact, the CAC model helps to 
obtain better estimates for some of these quantities and 
these estimates can then be used to generate better rate 
estimates. Using these estimates, CAC seems to produce 
KBB formation time scales which are not unrealistically 
long and, using our estimated parameters, are, in fact, 
on the order of lO"'* or 10"'^ per year or better. 

Another problem has be e n id entified by 
IChiang et al] ([2006) and iNoll et al.1 l|2007l ) with the L^s 
and mechanisms; assuming that the encounter time 
of the primaries in the Hill sphere tone ~ Rh/V then, 
for velocity dispersions with v > vh the time scale for 
binary formation goes as (n/uff)^^. On the other hand, 
if velocity dispersions are on the order of vh/3, i.e.. 



V < Vh, SiS assumed by IColdreich et al.l l|2002l ). then 
big bodies should coll apse into a vertically think disk 
(|Goldreich et al.ll2004 l with mutual orbit normals par- 
allel. This would imply that binary inclinations should 
be similar when, in fact they seem to be randomly 
distributed. 

In this section we compare CAC formation rates us- 
ing similar " back-of-the-enve l ope" a rguments to those 
provided by IColdreich et al.1 (|2002l ) and examine the 
fine-tuning problem. It is important to stress, however, 
that this discussion is not an attempt to compare the 
two models based on raw predictions of formation rates. 
The main aim is to deter mine if CAC lead s to ra tes sim- 
ilar to those obtained by IGoldreich et all (|2002D which, 
numerically, seem reasonable. In particular, based on 
our simulatio ns we conclude that sev eral assumptions in 
the model of IGoldreich et al] (|2002| ) may lead to over- 
estimates of binary formation rates. 



5.3 Binary formation rates 

Formation rates can be estimated in the L^s and L"^ 
models by making a few apparently reasonable assump- 
tions. Both processes (and also the CAC mechanism) in- 
volve the rate at which two big bodies enter their mutual 
Hill sphere. Assuming that the debris disk is infinitely 
thin, then the entry rate per big body c an be estimated 
as being th e product of three factors (|Goldreich et al.l 
|2002| . 120041 ): 

e « (^^^ X vug X Rh (19) 

where E is the surface mass density of big bodies of mass 
M. The first term is the number of big bodies per unit 
area. The second term is the velocity while the third 
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is an estimate of the range of impact parameters over 
which entry into the Hill sphere occurs. 

If it is assumed that bodies approach the Hill sphere 
with approximately the Hill velocity then eq. (|19|l be- 



e ; 



dR 



(20) 



where eq. (Ilip has been used. Here, d is the big body 
density and 6q is essentially a scaling parameter which 
relates the radius of a big body to the Hill sphere radius 
(|Goldreich et al.ll2004 ): i.e., Oq = R^iJ Rh ~ 10""; it 
turns out to be t he angular size of the Sun as seen from 
the Kuiper belt (|Goldreich et al.ll2004) The overall for- 
mation rate is the product of the Hill entry rate (O) 
and the probability (P) that an encounter with an in- 
truder leads to capture. Assuming prompt or ballistic 
transport through the Hill sphere leads to an average 
transient binary lifetime of tiij^ ~ Rh /Vh ~ 
The probability that a third large body will enter the 
quasi-bound binary Hill sphere during its lifetime is, 
therefore, 



dRhie 



1 

0% 



(21) 



The time scale for binary formation based on these two- 
body types of estimates becomes 



T2B 



pPQ 



dRhi 



2 a4 







(22) 



where p is the probability that, having entered the Hill 
sphere, a capture will occur. Setting p = 1 a nd us ing 
the paramete r s sugg ested bv lGoldreich et all (j2002l ) or 
IChiang et al.l (|2006l ) leads to a time scale on the order 
of 2 - 3 Myr for a quasi-bound binary to fuse. 

There are several difficulties with the approach 
just described, particularly in regard to computing the 
binary lifetime inside the Hill sphere. Binaries which 
transit the Hill sphere promptly, i.e., on a time scale 
« Rh /Vh ~ 50 years are rarely captured. Figure [2] 
shows that different trajectories having the same asymp- 
totic velocity can spend very different times inside the 
Hill sphere depending on the particulars of their orbital 
elements. Further, small asymptotic velocities need not 
imply slow transit of the Hill sphere (see Fig. [2] - the 
actual velocity after entering the Hill sphere can be 
substantially higher than the asymptotic value. While 
"rapid transit" binaries do indeed contribute to the 
overall entry rate into the Hill sphere their lifetimes are 
not relevant for computing the encounter rate of quasi- 
bound binaries with intruders. This is because chaos can 
induce a substantial time delay in the transit time of a 
fraction of the entering binaries. The lifetimes of these 
"resonant" states can be several orders of magnitude 
larger than is implied by prompt transport as Fig. [2] 
shows. Because their formation is a non-linear process 



these lifetimes cannot be predicted in any simple way. 
The situation is somewhat similar to a chemical reaction 
in which resonances - essentially the time-delay due to 
the formation and eventual decay of a quasi-bound state 
- can have a profound effect on the reaction cross-se ction 
(jLevine fc Bernsteinlll987l : lAlthorpe fc Clarvll20oi ). 

The possibility of very long-living resonant quasi- 
bound binaries suggests that an alternative - though 
still "back-of-the-en velope" style (|Goldreich et al.ll2002l : 
IChiang et al.l |2006| ) - treatment of binary formation 
rates might be useful. 

Assuming that the overall entry rate is roughly 
given by eq. p0|l we need to estimate the fraction of 
entrants which get caught up in chaotic layers as well as 
their lifetimes therein - i.e., their Hill lifetime. Suppose 
that the probability of a trajectory having a lifetime 
equal to, or greater than, some threshold Hill lifetime 
To is P>. Then the combined probability that a third 
large body will enter the quasi-bound binary Hill sphere 
and capture the binary permanently during its lifetime 
T > To becomes 



dRi 



X Py X To X p X Vl 



(23) 



where Ri is the radius of the intruder and vi is the in- 
truder velocity in units of Vh', we take vj ^ 2vh - see 
Fig. [1 This is the estimate of iGoldreich et al.1 l|2002h 
multiplied by a correction factor due to the delay inside 
the Hill sphere due to CAC. Although, in fact, p < 1, in 
or der to make a dir e ct com pari son with estimates m ade 
by IGoldreich et al. I (|2002t ) and IChiang et all (I2006D we 
also set p = 1 and Ri = -Rbig- This likely underesti- 
mates the actual capture time scales in the two-body 
simulations. 

The cumulative probability P> can be calculated 
directly from our simulations and thus the increase or 
decrease in the timescale for binary formation as com- 
pared to eq. (|22p can be estimated as a function of Tq. 
Figure [7] shows the ratio of the thusly computed CAC 
time scale (Tcac) and the corresponding two-body esti- 
mate based on eq. (|22p . i.e., T2B. We conclude that the 
time scales for binary formation in the CAC model are 
on the order of a few Myr. 

So far the discussion has neglected the small- 
intruder effect apparent in Fig. [5] If we assume that 
the intruder has a mass ~ 2% of the total binary mass 
then Ri ~ ii/3 and, for long-lived transient binaries we 
estimate from the simulations that p ^ 0.3 which, to- 
gether, decreases the time scale for capture by a further 
factor of ~ 3. Of course, decreasing p in eq. (18) seems 
to increase the CAC tim e scale compared to estimates 
bv lGoldreich et al.1 (|2002l ). However, it should be noted 
that in their work p was set to unity - i.e., the assump- 
tion was made that all intruder scattering events would 
lead to capture. Our results indicate that p < 1 and, 
for binaries which transit the Hill sphere rapidly - i.e.. 
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Figure 7. Ratio of two different estimates of the time scale for binary formation: Tqac is the CAC time scale and is estimated 
usi ng eq. I I23I I: T2B ~ 2 — 3 Myr is the corresponding two-body estimate from eq. ((22} ta oldreich et al.ll2002l : IChiang et al.l 
To is the threshold Hill lifetime (logarithmic scale) - see discussion surrounding eq. II23II . 




those which do not become entangled in chaotic layers 
- p « 1. Thus, in principle p should also be reduced 
when computing the two-body rates. 

Nevertheless, the CAC model, even though it per- 
mits supra-Uff velocity dispersions, and relies on the 
formation of relatively rare long-lived transient binaries, 
leads to capture rates which seem reasonable, i.e., on the 
same order-of-ma gnitude time scale as the growth of iso- 
lated big bodies iColdreich et al.ll200^ : lRafikov|[2003dl : 
iGoldreich et al.ll2004 ). 



5.3.1 Role of the Hill velocity 

Several observations on the significance of the Hill ve- 
locity in these simulations (and rate estimates) can be 
made. 

(i) The Hill velocity is defined to be the orbital ve- 
locity around a large central object at the Hill sphere 
radius. The assumption in analytical o r semi-analytical 
estim ates of binary formation rates (jColdreich et al.l 
I2OO2I ) is that if binary components are separated by 
« Rh then they orbit each other with the same pe- 
riod as th e binary center of mas s orbits the Sun, i.e ., 
l|Goldreich et al.l |2004| : IChiang et all |2006| ). 
However, this will rarely be the case, especially if the 
binary partners are following a chaotic orbit: galleries 
of orbits which inhabit outer regions of the Hill sphere 
MurisonI igSSl) together with their surfaces of section 
Winter fc Murravl1l994al lbh are available and these in- 



dicate a broad range of possible orbits, few of which 
seem likely to satisfy the Hill velocity criterion. In fact, 
the stability radius for two-body mutual moti on inside 
the Hill sphere is considerably less than Rh (llnnanen 
19791: lTothlll999l : [Nicholson fc Gladmanll2006l : ISzebehelv 



Therefore, it is unclear to what extent the Hill veloc- 
ity - an essentially two-body quantity - is useful in this 
situation. 

(ii) In view of the previous observation, and from the 
simulations, even if the velocity at infinity is ~ vh then 
this does not imply that the entry speed of the primaries 
- or the intruder - into the Hill sphere will be comparable 
to Vh. In addition, the instantaneous relative velocity 
of the quasi-bound primaries can be significantly larger 
than Vh. 

(iii) After the particles enter the Hill sphere what 
happens next is inherently unpredictable except by ex- 
plicit computation. Primaries which become entangled 
in chaotic layers may have very long lifetimes whereas 
others may transit promptly and it is difficult, in a 
chaotic system, to distinguish the two classes a pri- 
ori. It seems likely that more sophisticated phase space 
transport the ories will be need e d to e s timat e rates semi- 



analyticallv (IJaffe et all Il999l . |2000| . I2OO2I : lUzer et all 
I2OO2I : IWaalkens et al.ll2005bD . In principle, these rates 
will be related to the relative volumes of chaot ic and 
scattering parts of phase space (|jaffe et aLlbOOOD . 



1967); although stable orbits can be found at very large 
fractions of Rh their motion is highly non-Keplerian. 



5.4 The fine-tuning problem 

IGoldreich erall l|2002h have estimated that the veloci- 
ties of ~ 100 km sized bodies were on the order of vh /3 
in t he early Ku i per-b elt. However. IChiang et all l|2006l ) 
and I Noll et al] (|2007l ) note that if the velocity disper- 
sion, V, of big bodies is less than vh then these bodies 
will collapse into an essentially two-dimensional disk due 
to runaway cooling which is caused by the dynamical 
friction exerted by the sea of planetesimals. If the pri- 
mordial binary partners inhabited a vertically thin disk 
then TNB mutual orbit inclinations ought to be similar, 
i.e., mutual orbit normals should be approximately par- 
allel. In fact, TNB inclinations appear to be randomly 
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distributed ||No1]| l2003l : IChiang et all 120061 ) which im- 
plies an isotropic primordial velocity distribution. This 
suggests that the big bodies did n ot all originally lie in 
the same plane (|No11 et all [20071 ) and, therefore, that 
binaries formed when v > Vh- However, assuming that 
V > Vh has an unex pected snag. 

iNoll et al.1 (|2007l ) have demonstrated that rate es- 
timates based on two-body dynamics predict an ex- 
tremely sharp decrease in capture probabilit y for veloc- 
ity di spersions which exceed vh even slightly (|No11 et al.l 
|2007| ) - according to their work, capture time scales go 
as (v/vhY^ which essentially preclu des binary forma- 
tion unless v ^ Vh. INoII et al] |20o3) obtain this result 
by making the following assumptions: 

(a) For i; > vh the actual "sphere-of-influence" of 
the primaries is not Rh but Ri where 

Ri^RhX C^f (24) 
v 

(b) For V > Vh the average encounter time (in 
years) inside the sphere-of-influence is tcnc ^ Ri /V 

As velocities increase the sphere-of-influence shrinks 
as does t^nc', powers of v/vh proliferate and capture time 
scales ~ {v/vhY^ result. Obviously, according to this es- 
timate, even a small increase of v above vh will lead to 
dramatically long time scales for bina ry formation. Thi s 
constitutes the flne-tuning problem (INoU et al.l |2007| ): 
observed inclinations suggest that v > vh but rate com- 
putations indicate that v vh- Together, this implies 
that V = Vh wh ich is an example of fine-tuning. As 
INoU et all (|2007l ) ask "Why should v ^ vh during bi- 
nary formation?" There seems to be no good reason. 

The CAC model does not suffer from a flne-tuning 
problem which, at its heart, springs from the assumption 
that the encounter lifetime goes as Rh/V. This is for 
two reasons: 

(a) Figure [5] shows that, as predicted bv lNoU et al.l 
l|2007l ). the "sphere of influence" is reduced in size as 
compared to the Hill sphere radius. Intruders which 
penetrate deeper have a higher probability of lead- 
ing to capture. Nevertheless, the steepest drop-off in 
density in Fig. [S] occurs just outside Rh which sug- 
gests that setting Ri = Rh x (vh/v)^ (when v > vh) 
||No11 et allbOO'il ) is probably too severe of an approx- 
imation. Figure |5] shows the distribution of intruder 
asymptotic velocities leading to capture. These veloc- 
ities, typically, are comparable to or exceed vh- How- 
ever, small intruders in the primordial Belt were likely 
following more elliptical orbits and therefore their ve- 
locities may have substantially exceeded vh- Simula- 
tions in the elliptical Hill four-body problem will be 
needed to model such encounters. 

(b) In CAC the encounter lifetime, tone need not 
bear much relation to ?; - see Fig. [5] - and, in fact. 



the encounter lifetime inside the Hill sphere can be 
significantly extended by capture into chaotic layers. 
Thus, CAC is consistent with randomly distributed 
binary inclinations because it allows significant cap- 
ture probabilities when v > vh, i.e., it does not suffer 
from a flne-tuning problem. 

From these considerations we conclude that the CAC 
mechanism is, in principle, consistent with observed 
TNB mutual orbit inclination distributions. Because the 
Hill lifetime can be quite different from Rh/v it seems 
likely that the origin of the flne-tuning problem lies 
in the essentially two-body estimates which constitute 
its basis - it appear s that the fine-tuning problem was 
raised bv lNoll et al.l l|2007t ) precisely to make this point. 

5.5 Mass selection: Comparison with L'^s and 

The prediction of a preference for similarly sized binary 
partners in the CAC model departs from both the L^s 
dynamical f riction formation model and the L'' model 
(jColdreich et al. 200 2). 

L^s: For the small (relative to the Sun) masses of 
the TNB partners the Hill approximation is extremely 
accurate. In fact, as noted, the individual masses of the 
two small bodies - and the total mass - can be scaled 
out o f the equations of mot i on entirely (IPetit fc HenonI 
Il986l : iHenon fc Peti3 1 19861 : iMurrav fc DermottI llQQgT 
Therefore, simple dynamical friction cannot select for 
any particular mass ratio; in order for dynamical friction 
to be effective it suffices that the total mass of the poten- 
tial binary be much larger than the individual masses 
of planetesimals which make up the sea of small bodies 
exerting dynamical friction. This can be seen by not- 
ing that the semi-major axis which undergoes steady 
reduction through dy namical friction in the model of 
IColdreich et al.l (|2002l ) is actually the semi-major axis 
of the mutual orbit in the Hill approximation . Thus, 
while the L^s model of IColdreich et al.l (|2002l ) seems 
to require large primaries, it is only necessary that the 
total binary mass is large as compared to planetesimal 
masses. Therefore, this model does not make predictions 
about the mass ratio of the binary partners. 

L'': Although our calculations have demonstrated 
that three-body scattering inside the Hill sphere selects 
for roughly same-sized binaries it is unlikely that this es- 
sentially non-linear effect is predictable from two-body 
considerations. This is because chaos plays a double 
role: on the one hand it acts as the glue to hold together 
transient binaries for long periods of time, thereby en- 
hancing their effective capture cross sections. On the 
other hand, if the intruder itself happens to become en- 
tangled in a chaotic zone during its transit through the 
Hill sphere - an event that is more likely if the primaries 
have unequal masses - then the binary tends to be either 
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Figure 8. Intruder velocity distribution leading to capture. Units scaled by the Hill velocity vu. 



not s tabilized or actively destabilized (|Astakhov et alj 
I2OO5I ). Thus, the CAC mechanism discriminates not 
only between "big" and "small" bodies but, actually, 
also between large and small "big" bodies. 



6 BINARY HARDENING 

The question naturally arises as to under what condi- 
tions the intruder, on average, hardens - i.e., increases 
the binding energy - of the binary. We do not address 
this issue in detail here; for now, we note that vari- 
ous rules of thumb have been developed to understand 
whe n a binary hardens or softens in star-binary scatter- 
ing (Hut & Bahcall"l983';'Hut"l983';'Hiiyil983a"b','l98i, 
1990; HcEKic & Hut 2003; Aarscth 2003; Hcggic & Hu3 
19931: IHeggie et all 119991 For smaU intruders, iHiiy 
(1990) has argued for a "fast-slow" rule in which the 
binding energy of a binary increases (i.e., it hardens) if 
the velocity of the intruder at infinity is less than the or- 
bital velocity of the binary, and vice-versa. While large 
bodies in the disk generally are expected to have smaller 
velocities than small bodies, in the Hill sphere the bi- 
nary partners are initially following a chaotic trajec- 
tory influenced by solar tides. We flnd that quasi-binary 
partner relative velocities are, on average, higher than 
the asymptotic velocities of incoming intruders. This is 
due primarily to solar tides (which induce the chaos 
in the first place) . This is illustrated in Fig. |9] which 
shows the coordinate space projection of a very long 
lived trajectory. Frames (b) and (c) illustrate that the 
relative orbital frequency is substantially higher than 
the heliocentric frequency and, correspondingly, instan- 
taneous (and average) relative velocities can be substan- 
tially higher than asymptotic velocities. Therefore, on 
average, intruder scattering is expected to harden bina- 



ries u nless resonant triples (|Heggie fc Hutl2003l : lAarsethl 
I2OO3I ) form in which case the binary tends to be desta- 
bilized. This suggests that sequential scattering by rel- 
atively small intruders can be an efficient means of or- 
bit reduction, consiste nt with our previous simulations 
(|Astakhov et al.ll2005l ). 



7 LIMITATIONS OF THE CALCULATIONS 

The calculations reported here have not primarily been 
an attempt to calculate binary formation rates but, in- 
stead, are an investigation into the consequences of the 
CAC model. Several limitations of the current set of 
calculations are apparent. 

(i) The calculations assume that a quasi-bound bi- 
nary has already formed prior to the intruder being 
launched from inflnity. Because the mass and orbital 
element distributions at inflnity are poorly character- 
ized then certain attributes of the intruder were selected 
somewhat arbitrarily, e.g., the initial distance of the in- 
truder from the binary barycenter. While we do not ex- 
pect the qualitative conclusions of the calculations to be 
affected, the probability of binary formation will depend 
on these details. Ideally, initial masses and orbital ele- 
ments would be selected from a variety of possible mass 
and velocity distributions and the results compared with 
each other and with observations. Work to that end is 
in progress. 

(ii) A hierarchical separation of time scales is as- 
sumed in which the binary forms in the three-body Hill 
approximation and later intruder scattering occurs in 
the four-body Hill approximation. In principle, quasi- 
binary formation and intruder scattering could hap- 
pen essentially simultaneously. However, we expect that 
such encounters will be relatively infrequent. 
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Figure 9. An example of a proto-binary trapped in a chaotic layer. Frame (a) shows the i;,T] projection of the orbit. Frame 
(b) shows the radius p as a function of time scaled by the orbital period. All distances are scaled by the Hill sphere radius R[j. 
Frame (c) is a blow-up of part of Frame (b). 
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(iii) In these calculations exchange reactions were ig- 
nored and the intruder mass was not allowed to exceed 
the larger of the two primaries. Relaxing these assump- 
tions may lead to modifications of the observed distribu- 
tions; preliminary investigations indicate that exchange 
only tends to accentuate the preference for large mass- 
ratio binaries. 



8 CONCLUSIONS 

The model p r esent ed in our previous article 
(|Astakhov et al.l Hooi), and studied in more detail 
here, relies on the ability of chaos to produce very 
long-lived, quasi-bound binaries which can be stabilized 
by gravitational scattering with an intruder. Roughly 
equal-sized binary partners are predicted to form 



which translates into 



< 1. In addition the model 



suggests that binary formation occurs most readily 
when the intruder mass is a small percentage of the 
total binary mass. 

Several extensions of the model are possible. For ex- 
ample, it would be interesting to examine the effect of 
using particular distributions fr om which to pick initial 
masses and orbital elements (IGreenzweig fc Lissauej 
1990. 1 9 9i iKenvon fc Luul Il998l . Il999l : iKenvon 2002; 
Kenvon fc Bromlevll2004l ') : however, our numerical sim- 



ulations indicate that a trajectory started randomly in 
the feeding zone has ~ 0.8% chance of actually en- 
tering the Hill sphere - and a smaller probability of 
subsequently becoming caught up in a chaotic layer. 
This makes numerical simulations very inefficient in that 
most integrations result in no possibility of capture. 
That is, it is difficult, in practice, to select initial con- 
ditions directly in the asymptotic region because most 
do not lead to interpenetration of the Hill sphere. 

To address these issues and to increase the effi- 
ciency of the scattering calculations we are in the pro- 
cess of developin g a method which u ses random forests 
of decision trees l|Breiman et al ]|1984| ) as a Boolean clas- 
sifier to discriminate between instances - sets of initial 
conditions at infinity - which lead to primary encoun- 
ters inside the Hill sphere (positives) and those which 
do not (negatives). 

To understand the capture dynamics further and 
to predict formation rates it m ay be fruitful to employ 
phas e -space transport theory (|jaffe et al.l 1 19991 . |2000| . 
l2002l : IUzer et al.ll2002h to study both the initial tempo- 
rary trapping of the primaries and their subsequent sta- 
bilization or destabilization through intruder scattering. 
Essentially the procedure is as follows: (i) identify the 
critical saddle-points (in the Hill problem these are the 
Lagrange points L\ and L2), (ii) make a Birkhoff nor- 
mal form expansion at the saddle points, (iii) construct 
the normally invariant hyperbolic man ifolds and the 
transition state(s) from the normal form l|Wiggins et al.l 



I1994I ). The transition state is a dividing surface in phase 
space that acts a surface of no return. This method 
has already been applied to the three-body Hill prob- 
lem, although o nly over a quite limited energy range 
(jWaalkens et al.i i2005b£iai) • 

New observations will provide rigorous tests of 
the CAC scenario. Observationally, few of the already 
known TNBs are fully characterized. In particular the 
true inclination of the mutual orbital plane generally 
remains unknown. Because of the small parallax, am- 
biguities in the orbit inclination of TNBs will per- 
sist unt il observations are performed at well separated 
epochs (iHe stroffer fc Vach icr 2005). Knowledge of true 
inclinations is especially important for predicting stel- 
lar occultations and mutual p henomena from which bi- 
nary densities can be derived (|Brown fc Webster! [1997I : 
iPravec et al] 1 19981 : ICooravl |2003| ). Orbital inclinations 
are of interest here for two reasons: (i) the CAC for- 
mation model predicts asymmetric distributions of pro- 
grade vs. retrograde orbits with respect to the invariant 
plane; (ii) the numerical results presented here and ac- 
tual observations of inclination distributions are consis- 
tent with TNB production during the oligarchic growth 
phase which comes directly aft er a shor t er, run away ac- 
cretion regime in the model of lRafikovl l|2003al '). 

In summary, trans-Neptunian binaries present mul- 
tiple conundrums which include their unexpectedly high 
proportion in the total population of TNOs; their large, 
mid-range eccentricity, orbits; their mutual orbital in- 
clination distributions; and their apparent propensity 
to have roughly equal-sized partners. Resolving these 
questions will undoubtedly shed light on the forma- 
tion of the trans-Neptunian part of the Solar Sys- 
tem, with i mplications for plane t ary formation theories 
in general (lYoudin fc ShulboOSl : IColdreich et"aLl 12004 
iBasri fc Brownll2006l ). The results presented here sug- 
gest that chaos might have played an important - and 
constructive - role in the formation of these remote and 
puzzling objects. 
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APPENDIX A: KERNEL SMOOTHING 

Scatterplots are a useful way of visualizing data, e.g., the 
data in Fig. 2 could be represented simply as a scatter- 
plot of logj^Q(THiii) against velocity. However, when there 
exists a very high density of points in certain localized 
regions it may become difficult to dist i nguish local den- 
sity gradients llCentleman et al.|[2006l : ICarr et aLlll987l : 
lGentlemanll2007l l. In kernel smoothing the data is di- 
vided into two-dimensional bins and the bin coloured 
according to the local density using a false-colour rep- 
resentation. A kernel density estimator is used to es- 
timate the local density of points and to remove the 
dependence on the end points of the bins. The esti- 
mator centers a kernel function at each data point, 
e.g., a Gaussian. The kernel estimator essentially smears 
the contribution of each data point over a local re- 
gion which depends on the bandwidth and which, in 
turn, determines the relative contributions of adjacent 
points. Visual comparison of the raw scatterplot at 
high resolution and the kernel-smoothed plots shows 
that kernel smoothing enhances nuances in the den- 
sity variation that are hard to see in the scatterplot 
at the same resolution. Standard routines are available 
in the R statistical package for a n alysis and graphics 
iIR Development Core TeamI 12 007: Ha rrell et all l2006l : 
lNeuwirthl l2005l : I Gentleman et al.. ,2004 ) to make these 
plots. Because this is primarily a visualization proce- 
dure it is difficult to assign a quantitative colour scale; 
just as in a conventional scatterplot the density of points 
is best judged by eye. 

A simple analogue might be using a moving aver- 
age to smooth a time series in one-dimension and then 
making a histogram of the moving average rather than 
of the raw data. 
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